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5_i . Abstract 

We test the effectiveness of the multicanonical algorithm for the tertiary struc- 
ture prediction of peptides and proteins. As a simple example we study Met- 
enkephalin. The lowest-energy conformation obtained agrees with that determined 
by other methods such as Monte Carlo simulated annealing. But unlike to simulated 
annealing the relationship to the canonical ensemble remains exactly controlled. 
— |- . Thermodynamic quantities at various temperature can be calculated from one run. 

A protein or a peptide is a molecule that consists of a chain of N amino acid residues. 
There are 20 different amino acids known in nature. When N is large one calls the 
molecule a protein, otherwise a peptide. The prediction of tertiary structures of proteins, 
which determine their biological function, from their primary sequences remains one of 
the long-standing unsolved problems (for recent reviews, see, for example, Refs. 0). It is 
widely believed that this structure corresponds to the global minimum in the energy. So 
■ the problem amounts to finding the global minimum energy out of a huge number of local 

minima separated by high tunneling barriers. Within the presently available computer 
q ! resources, the traditional methods such as molecular dynamics and Monte Carlo simula- 

tions at relevant temperatures tend to get trapped in local minima. One of the methods 
which which seem to alleviate this multiple- minima problem is simulated annealing. 
However, a disadvantage of simulated annealing is that there is no established protocol 
for annealing and a certain number (which is not known a priori) of runs are necessary 
to evaluate the performance. Moreover, the relationship of the obtained conformations to 
the equilibrium canonical ensemble at a fixed temperature remains unclear. 



These problems may be overcome by the multicanonical algorithm which was recently 
proposed by Berg et al.M Originally developed to overcome the supercritical slowing 
down of first-order phase transitions, |4[] it has also been tested for systems with con- 
flicting constraints such as spin glasses. || H, |?J The latter systems suffer from a similar 
multiple-minima problem and it was claimed that the multicanonical algorithm outper- 
forms simulated annealing in these cases. |J 
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The idea of this method is based on performing Monte Carlo simulations in a multi- 
canonical ensemble 0, |] instead of the usual (canonical) Gibbs-ensemble. In the canoni- 
cal ensemble, configurations at an inverse temperature (3 = 1/RT are weighted with the 
Boltzmann factor Vb(E) = exp (^—(3E^j . The resulting probability distribution is given 
by 

P B (E) oc n(E)V B (E) , (1) 

where n(E) is the spectral density. In the multicanonical ensemble, |8| on the other 
hand, the probability distribution is defined in such a way that a configuration with any 
energy enters with equal probability: 

P mu (E) oc n(E)V m u(E) = const. (2) 

Then it follows that the multicanonical weight factor should have the form 

V mu (E) oc n-\E) . (3) 

In order to define a explicit form of this weight factor, we introduce two parameters a(E) 
and (3(E) as follows: [1 

V m u{E) = exp{-(/3 + (3(E))E - a(E)}. (4) 

For any fixed /3(E) and a(E) this leads to the canonical weight factor with the inverse 
temperature (3 = (3 + (3(E), hence the name "multicanonical". For a numerical simulation 
one needs estimators for the multicanonical parameters (3(E) and a(E). The iterative 
procedure by which one can get such estimators is described elsewhere. Once the mul- 
ticanonical parameters are determined, one multicanonical run is in principle enough to 
calculate all thermodynamic quantities by re- weighting. | T0|] Since in the multicanonical 



ensemble all energies enter with equal probability a simulation may overcome the barriers 
between local minima by connecting back to the high temperature states. In this way the 
global minimum can be explored. 

In the present work we apply the multicanonical algorithm to the problem of pro- 
tein folding, the tertiary structure prediction of peptides and proteins. The purpose of 
this work is primarily to test the effectiveness of the algorithm. For this reason we have 
studied one of the simplest peptide, Met-enkephalin. The lowest-energy conformation for 
the potential energy function ECEPP/2 Jll]| is known fl~2|1 and analyses with Monte Carlo 



simulated annealing with ECEPP/2 also exist. ||T3|, [14] 



Met-enkephalin has the amino- acid sequence Tyr-Gly-Gly-Phe-Met. For our simula- 
tions the backbone was terminated by a neutral NH2- group at the N-terminus and a 



neutral -COOH group at the C-terminus as in the previous works of Met-enkephalin. |T2 
The potential energy function that we used is given by the sum of the electrostatic term 
E es , the van der Waals energy E vdW , and hydrogen-bond term E hb for all pairs of atoms 
in the peptide together with the torsion term E tors for all torsion angles: 

Etot = E es + E v dw + Ehb + Etors (5) 

E„ = (6) 
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Figure 1: Definition of dihedral angles 
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Figure 2: Average energy < E > and specific heat C of Met-enkephalin as a function of 
temperature evaluated by multicanonical algorithms. The results of canonical simulations 
at fixed temperatures (50 K and 300 K) are also plotted (□). 
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is the distance between the atoms % and j, and ot\ is the torsion angle for the chemi- 
cal bond I. For a definition of these angles which represent the true degrees of freedom 
see Fig. |I]. The parameters Ay, Bij,Cij, Dy,Ui and nj) for the energy function were 
adopted from ECEPP/2,JTTJ]. The effect of surrounding atoms of water is neglected and 
the dielectric constant c is set equal to 2. The computer code KONF90,||rS| was modified 
to accommodate the multicanonical method. The peptide-bond dihedral angles uj were 
fixed at the value 180° for simplicity which leaves 19 angles 4>u ^% an d X% as independent 
variables. 



In Fig. 2a we show the average energy obtained by our method in a run with 10 5 



Table 1: Table I. Energy and dihedral angles of the lowest-energy conformations of Met- 
enkephalin obtained by multicanonical runs. Conformation A is the lowest-energy con- 
formation obtained by Monte Carlo simulated annealing (taken from Ref. 13). 



Conformation 


A 


1 


2 


3 


4 


5 


6 


E [ kcal/mol ] 


-11.9 


-11.9 


-12.0 


-12.0 


-12.1 


-12.0 


-11.9 


01 


98 


90 


91 


90 


97 


96 


98 




154 


153 


152 


154 


151 


153 


156 


02 


-161 


-160 


-157 


-161 


-158 


-161 


-163 


"02 


69 


72 


64 


71 


71 


68 


65 


03 


65 


64 


66 


63 


64 


64 


66 


■03 


-93 


-95 


-92 


-95 


-94 


-89 


-92 


04 


-85 


-82 


-80 


-77 


-83 


-85 


-80 


■04 


-27 


-26 


-29 


-32 


-30 


-31 


-29 


05 


-83 


-81 


-82 


-78 


-80 


-82 


-86 


^5 


142 


142 


138 


137 


145 


151 


147 


xl 


-179 


179 


-177 


179 


179 


-178 


-176 


xi 


-112 


-110 


-117 


-109 


-111 


-115 


-114 


xl 


149 


144 


146 


143 


149 


145 


142 


xl 


180 


-176 


178 


177 


180 


-178 


180 


xl 


73 


79 


81 


86 


79 


78 


78 


xl 


-65 


-64 


-67 


-67 


-66 


-67 


-66 


xl 


180 


-179 


180 


180 


-176 


180 


176 


xl 


179 


178 


179 


-179 


-179 


-178 


-178 


xi 


-55 


-66 


-59 


-62 


-61 


-60 


-57 



sweeps, as a function of temperature. The value ~ —12 kcal/mol at T = 50 K is very 
close to the global-minimum energy obtained by other methods. [jP2"| , |13| , H] In Fig. 2b we 
likewise present the "specific heat" (per residue), which is defined by 

C = f3 2 . (10) 

It has a peak around T = 300 K, which indicates that this temperature is important for 
peptide folding. 

During the production run the system reached the global-energy minimum region six 
times. The lowest-energy conformation within each visit is listed in Table I together 
with the global-minimum energy conformation (Conformation A in Table I) obtained by 
simulated annealing. |y| Conformations 1-6 are the results at Monte Carlo steps 20128, 



39521, 44462, 65412, 89413, and 95143. Hence, the system reached the lowest-energy re- 
gion in every 5000 to 20000 Monte Carlo steps. The energies are almost all equal, and the 
lowest-energy value in the present work (—12.1 kcal/mol) is slightly less than the previ- 
ous result (—11.9 kcal/mol) by simulated annealing. |13| Most of the dihedral angles of the 
six conformations also agree with the corresponding ones of Conformation A within « 5° . 



We have applied the recently developed multicanonical algorithm to the problem of 
predicting the peptide conformation. This method avoids getting trapped in a local min- 



imum of energy function by connecting back to high temperature states and enhances in 
this way the probability to find the global minimum. We have demonstrated the effective- 
ness of the algorithm by reproducing the lowest-energy conformation of Met-enkephalin. 
Furthermore, the multicanonical algorithm can yield various thermodynamic quantities 
as a function of temperature from only one production run. 
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